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ABSTRACT 

Two recent empirical developments in the study of extragalactic globular cluster (GC) populations are the 
color-magnitude relation of the blue GCs (the "blue tilt") and the nonlinearity of the dependence of optical 
GC colors on metallicity. The color-magnitude relation, interpreted as a mass-metallicity relation, is thought 
to be a consequence of self-enrichment. Nonlinear color-metallicity relations have been shown to produce 
bimodal color distributions from unimodal metallicity distributions. We simulate GC populations including 
both a mass-metallicity scaling relation and nonlinear color-metallicity relations motivated by theory and 
observations. Depending on the assumed range of metallicities and the width of the GC luminosity function 
(GCLF), we find that the simulated populations can have bimodal color distributions with a "blue tilt" similar to 
observations, even though the metallicity distribution appears unimodal. The models that produce these features 
have the relatively high mean GC metallicities and nearly equal blue and red peaks characteristic of giant 
elliptical galaxies. The blue tilt is less apparent in the models with metallicities typical of dwarf ellipticals; the 
narrower GCLF in these galaxies has an even bigger effect in reducing the significance of their color-magnitude 
slopes. We critically examine the evidence for nonlinearity versus bimodal metallicities as explanations for the 
characteristic double-peaked color histograms of giant ellipticals and conclude that the question remains open. 
We discuss the prospects for further theoretical and observational progress in constraining the models presented 
here and for uncovering the true metallicity distributions of extragalactic globular cluster systems. 
Subject headings: galaxies: elliptical and lenticular, cD — galaxies: fundamental parameters — galaxies: star 
clusters — globular clusters: general 



1. INTRODUCTION 

Massive early-type galaxies have globular cluster (GC) sys- 
tems that follow bimodal color distributions (e.g., Gebhardt 
& Kissler-Patig 1999; Larsen et al. 2001; Harris et al. 2006; 
Peng et al. 2006, hereafter P06). Before this was a commonly 
recognized feature of such systems. Ashman & Zepf (1992) 
discussed the possibility of GC color bimodality based on the 
simple idea that ellipticals are the remnants of dissipationally 
merged pairs of spirals (e.g., Toomre 1977). In this picture, a 
merger-formed population of red, relatively high-metallicity 
GCs would then appear distinct from the blue, metal-poor 
GCs of the progenitor spirals. 

The actual formation histories of elUpticals appear to have 
been more complex (see for example the review by Renzini 
2006). In particular, most of the stars in large ellipticals today 
are believed to have formed early, in a series of rapid star- 
bursts at redshifts z > 2, whereas the hierarchical assembly of 
these galaxies continued to z < 1 through "dry" (dissipation- 
less) mergers (e.g., De Lucia et al. 2006; Menci et al. 2008; 
Kodama et al. 2007; Thomas et al. 2005). Such hierarchi- 
cal formation models do not naturally produce bimodal GC 
metallicity distributions in most elliptical galaxies (Beasley 
et al. 2002; Kravtsov & Gnedin 2005). Instead, the stochastic 
nature of the formation would be expected to result in approx- 
imately Gaussian distributions. Thus, while some galaxies 
with quiet histories like the Milky Way's may have bimodal 
GC metallicity distributions if they formed their stars in just a 
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few discrete bursts, it is difficult to see how bimodality would 
be nearly universal in this context. 

An alternative model was put forth by Cote et al. (1998) 
in which ellipticals form with unimodal, relatively rich GC 
metallicity distributions and later acquire the blue GC compo- 
nent through dissipationless mergers with early-type dwarfs. 
This explanation appears more consistent with hierarchical 
formation models. However, the universality of the bimodal- 
ity, even at intermediate luminosities, and the asymmetric tail 
of red GCs in dwarf ellipticals, would again be surprising in 
this scenario. 

Models of GC color distributions have typically assumed 
linear conversions between metallicity and optical colors. The 
observed color bimodality thereby became equated with bi- 
modality in metallicity. In contrast, P06 showed that the em- 
pirical transformation from metallicity to the broad-baseline 
ig-z) color is distinctly nonlinear An earlier indication of 
such nonlinearity for the {C-Ti) color index (Harris & Har- 
ris 2002; Cohen et al. 2003) has also been confirmed (Lee 
et al. 2008). Further, Yoon et al. (2006) pointed out that 
nonlinear color-metallicity relations predicted by their sim- 
ple stellar population (SSP) models naturally produce bi- 
modal color histograms from unimodal metallicity distribu- 
tions. Such nonlinearities are commonly predicted by realistic 
SSP models that try to reproduce the structure in stellar color- 
magnitude diagrams as a function of metallicity (see Cantiello 
& Blakeslee 2007, hereafter CB07). This provides an alterna- 
tive, simple explanation for the near-universality of bimodal 
GC color distributions, although the idea remains controver- 
sial as we discuss in Sec.|4]below. 

It is not surprising that as observations and models improve, 
deviations from the simplest linear assumptions become im- 
portant. A related example is the dependence of galaxy sur- 
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face brightness fluctuation magnitudes on integrated color: 
linear fits sufficed at the level of accuracy of ground-based 
data, but the latest HST studies (Mei et al. 2007; Blakeslee 
et al. 2009) clearly show a wavy relation reminiscent of 
that between metallicity and color Galaxy scaling relations 
are another example where nonlinearity has become evident 
amid the continuity of structural properties in new samples 
of galaxies spanning many orders of magnitude (Cote et al. 
2008). Even the Cepheid period-luminosity relation has been 
found to exhibit significant nonlinearity (Ngeow et al. 2009). 
Although the errors in derived quantities (e.g., distance modu- 
lus or metallicity) resulting from such nonlinearities are often 
small for individual objects, the implications for the distribu- 
tions of these quantities derived for large populations can be 
profound. 

Another development has been the discovery that GCs in 
the blue peak of the color distribution tend to be redder at 
brighter magnitudes (Harris et al. 2006; Mieske et al. 2006; 
Strader et al. 2006). Dubbed "the blue tilt," the effect is 
subtle, but further evidence confirms it is real (Harris 2009; 
Peng et al. 2009). The likely explanation is that the most 
massive GCs were able to self-enrich to a small extent, thus 
producing a mass-metallicity relation (e.g., Dopita & Smith 
1986; Morgan & Lake 1989; Brown et al. 1991; Recchi & 
Danziger 2005; Strader & Smith 2008; BaiHn & Harris 2009). 
The tendency is much less obvious for the red-peak GCs, de- 
spite the otherwise lockstep behavior of the two peaks, which 
is strong evidence that they formed under similar conditions 
(P06; Larsen et al. 2001). However, the red peak is broader 
because of the nonlinear color-metallicity relation, and it usu- 
ally does not reach as high in optical luminosity, at least partly 
because of the decrease in luminosity-to-mass ratio with in- 
creasing metallicity. A weak "red tilt" is expected theoreti- 
cally (Bailin & Harris 2009), and there is evidence for this in 
some cases (Mieske et al. 2006; Lee et al. 2008). 

The purpose of this work is to bring together these two re- 
cent developments in the study of extragalactic GC systems: 
we explore the nature of the blue tilt in the context of uni- 
modal metallicity distributions with nonlinear dependence of 
color on metallicity. To do this, we assume a mass-metallicity 
relation consistent with the findings of recent studies. The 
following section describes our approach in more detail, and 
Sec. |3] presents the consequences for the observed blue tilt 
under various assumptions about the range of GC metallici- 
ties and luminosities. The implications of our results and the 
question of bimodality are considered in Sec.|4l The final sec- 
tion summarizes our conclusions. 

2. METHOD 

Our goal is to investigate the presence and significance of 
the "blue tilt" — the color-magnitude relation of the blue 
peak GCs — in the case of unimodal metallicity distributions 
following an empirically motivated mass-metallicity relation 
(MMR) and the nonlinear color-metallicity relations both ob- 
served and predicted by detailed SSP modeling. We also ex- 
amine the "tilt" of the red peak GCs in this scenario. For the 
sake of generality in this initial investigation, we try to keep 
our assumptions as simple as possible while still being realis- 
tic enough to yield useful insights. 

2.1. Color-Metallicity Relations 

P06 presented an empirical [Fe/H]-(g-z) relation based on 
their photometry and spectroscopic metallicities for a com- 
bination of 95 Milky Way, M49 (NGC4472), and M87 GCs 
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Fig. 1 . — Empirical (g—z) and model (B-I) color-metallicity relations. The 
circles show data for Milky Way and Virgo cluster GCs (Peng et al. 2006), 
and the long-dashed curve is a fourth-order polynomial fit to the data. The 
diamonds show predictions for (B-I) from the 13 Gyr Teramo SPoT models 
used in this study. The solid curve is a fifth-order polynomial fit, while the 
short-dashed curve is a cubic spline interpolation/extrapolation of the models. 
As commonly used for these indices, the (B-I) colors are on the Vega-based 
Johnson-Cousins system while (g-z) is given in AB magnitudes. The con- 
version from the AB to Vega system for (g—z) is +0.64 mag; thus, these color 
indices are quite similar, but the data have a steeper slope than the models at 
the red end. 



(see their Fig. 11). The Milky Way GC data are part of a larger 
project to produce a homogeneous database of integrated col- 
ors for the entire Galactic GC system in the SDSS photometric 
bandpasses (West et al. 2006; M. West et al., in preparation). 
The inclusion of the Virgo galaxy GC metallicities from Co- 
hen et al. (1998, 2003) was necessary to characterize the re- 
lation at high metallicity. This adds uncertainty because the 
conversion from the Lick indices to the Zinn & West (1984) 
metallicity scale is not well calibrated above solar metallicity 
and involves some extrapolation (see Fig. 14 of Cohen et al. 
1998). However, the extrapolation is strictly linear and ex- 
tends < 0.5 dex at the metal-rich end, so it could not introduce 
a nonlinearity in the [Fe/H]-(^-z) relation. Intrinsic trends in 
age and/or alpha-enhancement with metallicity could poten- 
tially introduce biases in the metallicity scale (see the discus- 
sion in P06), but such trends would not change our results if 
they involve a simple scaling factor. Since we wish to avoid 
ad hoc assumptions about model ages and abundance ratios 
as far as possible, we prefer to adopt empirical calibrations 
in our analysis, and this remains the best empirically-based 
color-metallicity relation available. 

For simplicity, P06 used a broken linear fit to the curved 
[Fe/H]-(^-z) relation, excluding the highest metallicity GCs 
which would have required additional curvature. To better 
characterize the full data set presented in P06, we fit a quar- 
tic polynomial using robust orthogonal regression (GaussFit; 
Jefferys et al. 1988) and find the following relation 

[Fe/H] = -33.74 + 19M(g-z)-66.69(g-zf 

+ 20.66(g-zf-im(g-z)\ (1) 

which we plot in Figure [T] along with the data points from 
P06. Note that the empirical (g-z) colors are calibrated to 
the AB system, and the conversion from AB to Vega-based 
magnitudes for (g-z) is -1-0.64 mag. 

Other color indices used in recent studies of GC color dis- 
tributions and MMRs include (B-I) (e.g., Harris et al. 2006; 
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Fig. 2. — Calculated (g—z) color distributions based on the empirical poly- 
nomial fit in Fig. [T] for a series of unimodal Gaussian [Fe/H] distributions 
having a common dispersion of 0.5 dex and means ranging from -1.2 to 
-0.2 dex in steps of 0.1 dex. The normahzation of the curves is arbitrary; 
the color ranges from violet/blue at low metallicity to red/magenta at high 
metallicity. Although the mean color increases with metallicity and the peaks 
shift redward, the position of the trough remains fairly constant near the point 
where the derivative of the metallicity with respect to color has a minimum. 

Harris 2009), (V-I) (Peng et al. 2009), (B-R) (Spitler et al. 
2006; Cantiello et al. 2007), (g-i) (Wehner et al. 2008), and 
(C-Ti) (Forte et al. 2007; Lee et al. 2008). The empiri- 
cal metallicity dependencies of these indices have not been 
characterized in the same detail as (g-z)', usually the con- 
versions to metallicity are based on simple linear approxima- 
tions. However, CB07 present color-metallicity predictions in 
the UBVRIJHK bandpasses from the Teramo "SPoT" mod- 
el^ (Raimondo et al. 2005) and make comparisons to other 
SSP models in the context of color-metallicity relations. Fig- 
ure [U also shows the color-metallicity predictions from the 
13 Gyr SPoT models for (B-I), the color index which first 
revealed the presence of the MMR in the GCs of giant ellip- 
ticals (Harris et al. 2006), along with simple polynomial and 
cubic spline fits to the models. (The age of 13 Gyr is the 
same as for the models used by Yoon et al. [2006] and CB07.) 
The plotted (g-z) and (B-I) relations provide empirical and 
theoretical examples of nonlinear color-metallicity relations 
for use in this study. See Yoon et al. (2006) for a theoretical 
(g-z)-metallicity relation and CB07 for predicted relations in 
other bandpasses. 

2.2. Mass-Metallicity Scaling 

There is evidence that the MMR is already present near the 
peak of the GC luminosity function (GCLF), although it is 
more prominent a couple magnitudes brighter than this (e.g., 
Harris et al. 2006; Mieske et al. 2006; Peng et al. 2009). Using 
linear color-metallicity relations for blue GCs, these studies 
derive scaling relations between GC luminosity L and metal- 
licity Z consistent with Z ^ L*' for massive elliptical galaxies 
(e.g., Cockcroft et al. 2009; Strader & Smith 2008), or a slope 
of -0.2 dex mag"' for [Fe/H] (we confine the discussion here 
to solar-scaled models). This coefficient is not observationally 
well constrained, as the galaxy-to-galaxy scatter is large and 
the applied color-metallicity conversions are approximate. 

We implement this MMR scaling by assuming that the Z 
distribution function (ZDF) is a Gaussian in the logarithm 
with cr[Fe/H] = 0.5 dex (e.g., Yoon et al. 2006; CB07) at each 

' http://www.oao-teramo.inaf.it/SPoT 
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Fig. 3. — (B-I) color distributions given by the polynomial (left) and spline 
(right) fits to the model points in Fig. [T] for a series of unimodal Gaussian 
[Fe/H] distributions having a common dispersion of 0.5 dex and means rang- 
ing from -1.2 to —0.2 dex in steps of 0.1 dex. The normalization of the 
curves is arbitrary; the color ranges from violet/blue at low metallicity to 
red/magenta at high metallicity. Although the polynomial and spline fits to 
the model points in Fig. [T] are very similar, the derived color distributions 
depend on their derivatives and thus show marked differences. For real pop- 
ulations, observational errors would largely obscure such differences. 

magnitude, but with the mean [Fe/H] increasing by 0.2 dex 
per magnitude brighter than the GCLF peak, up to some limit 
(see Sec. 12.31 ). Figure|2]shows the predicted (g-z) color distri- 
butions for Gaussian ZDFs with o-fpe/Hj = 0-5 dex and means 
varying from -1.2 to -0.2 dex. Consistent with the obser- 
vations (e.g., P06), even at low mean metallicities there is a 
substantial "red tail" to the color distribution; moreover, the 
distribution always appears distinctly bimodal at the highest 
metallicities appropriate to GC systems, i.e., the blue peak is 
always present. 

Figure [3] presents the corresponding plots for (B-I) based 
on the polynomial and spline fits to the SPoT model predic- 
tions. We show both sets of curves to illustrate the fact that 
since the color distributions depend on the derivative of the 
nonlinear relation between metallicity and color, the interpo- 
lation method can affect the details of the derived distribu- 
tions. In practice, observational errors would obscure most of 
these details. For the analysis below, we concentrate on the 
results derived using the polynomial interpolation to be more 
consistent with the approach used for the empirical (g-z) rela- 
tion (for which spline interpolation is not a reasonable option). 
However, we have carried out the analysis using both meth- 
ods, and the results are qualitatively the same. In this case, 
the uncertainties and limitations of the model predictions are 
greater than the differences due to the interpolation scheme. 

The GCs in the simulations follow a standard Gaussian 
GCLF at any given metallicity. Since the MMR is most 
prominent in giant ellipticals, we use a GCLF width ctlf =1.4 
mag, characteristic of the most luminous galaxies (e.g., Har- 
ris 1991; Blakeslee et al. 1997; Harris et al. 2009), except in 
the discussion of the MMR for dwarf ellipticals in Sec. 13.21 
below. We apply the model mass-to-light values to obtain the 
correct relative luminosities at different metallicities, so that 
the overall population follows a single Gaussian distribution 
in mass (e.g.. Ashman et al. 1995). 

2.3. Monte Carlo Simulations 

The main adjustable parameters in this analysis are the 
mean [Fe/H] at the peak of the GCLF, which we label 
([Fe/H])o, and the maximum mean metallicity that we allow 
for the brightest GCs, which we label ([Fe/H])niax- In other 
words, we impose a fixed mean of ([Fe/H]) max for the ZDF of 
the GCs above some luminosity limit, rather than continuing 
the scaling relation indeterminately. The mean [Fe/H] of the 
GC system could be derived from ([Fe/H])o, ([Fe/H])niax, and 
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the GCLF. Empirically, the mean GC metallicity increases 
with galaxy luminosity (Cote et al. 1998; P06). We do not ex- 
plicitly include galaxy properties in these simulations, but we 
consider a large enough range in ([Fe/H]}() and ([Fe/H])n,ax 
to make useful comparisons to galaxies with a range of lumi- 
nosities. 

Given the above set of assumptions, we generate random re- 
alizations of GC systems and examine their color distributions 
and color-magnitude relations. We simulate realistic observa- 
tional scatter by adding random magnitude errors that scale 
exponentially as ct,„ ~ lO""'^™, and we set the normalization 
such that the color error at the GCLF peak is 0.1 mag. We 
also include an additional 0.01 mag scatter at all magnitudes, 
which could represent errors in aperture corrections or intrin- 
sic scatter in the GC colors at fixed metallicity. We experi- 
mented with larger values of these errors, but we note that the 
blue tilt is only detectable with the best available photometric 
data and analyses, and we do not wish to obscure the effect 
we are trying to characterize. 

2.4. Color-Magnitude Fitting Procedure 

To characterize the blue (and red) tilts in our simulated sam- 
ples, we basically follow the procedure described in Mieske 
et al. (2006), which is similar to that used in other studies 
(e.g., Peng et al. 2009). Specifically, we bin the GCs by mag- 
nitude, using a minimum bin size of 0.2 mag but allowing the 
size in magnitudes to increase in order to contain at least 200 
GCs per bin. We then use the heteroscedastic mode of the 
KMM mixture modeling algorithm (Ashman et al. 1994) to 
determine the peaks of the color distributions for the different 
bins. We perform linear fits to determine how the peak posi- 
tions vary as a function of magnitude. The fits are done using 
all bins that are more than 1 ctlf brighter than GCLF turnover 
magnitude; for giant ellipticals this is close to the "million so- 
lar mass" threshold (Harris 2009) where the tilt becomes no- 
ticeable in observational data. The following section presents 
the results for the color-magnitude slopes for the blue and red 
peaks in these simulations. 

3. RESULTS 
3.1. Blue vs Red Tilts for Baseline Models 

Figure |4] plots the slope of the blue peak in the ig-z) color 
distribution determined by KMM as a function of magnitude 
versus the corresponding slope for the red peak. Each panel 
shows results for 101 random realizations of a GC popula- 
tions with 10,000 members brighter than the GCLF peak. 
The twelve different panels are for different combinations of 
([Fe/H])o and ([Fe/H])max, defined in Sec. 12.31 their values 
are shown in each panel. In particular, the mean metallicity 
at the GCLF peak ([Fe/H])o is varied from -1.0 dex (one- 
tenth solar) to -0.5 dex (one-third solar), but only two val- 
ues, -0.3 (half solar) and 0.0 (solar), are shown for the max- 
imum mean metallicity parameter ([Fe/H])n,ax. Lower val- 
ues of ([Fe/H])n,ax decrease the mean metallicity range of the 
MMR to factors of only ^2 and/or restrict the simulated pop- 
ulations to lower mean metallicities and bluer mean colors 
than those found in the giant ellipticals where the "blue tilt" is 
commonly observed. Setting ([Fe/H])max to values above so- 
lar metallicity (or removing it altogether) gives results similar 
to the case with ([Fe/H])inax = 0.0 dex, already a high metal- 
licity for GCs. Thus these two values of ([Fe/H])max prove 
adequate for our purposes. 

The points in Figure |4] with negative slope values are the 
ones that show a reddening of the peak colors with increasing 
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Fig. 4. — Blue versus red tilts in simulated GC ig—z) color distributions. 
Each point represents the slopes of the blue and red peak colors determined 
by KMM fits as a function of magnitude within a single realization of a popu- 
lation of 10,000 GCs. See text for details on KMM fitting. Each panel shows 
the slopes for 101 different realizations (squares), all with the same model 
parameters. The large encircled crosses mark the average values. From panel 
to panel, we vary two parameters: ([Fe/H])o, which is the mean metallicity 
for GCs at the mean luminosity of the GCLF, and ([Fe/H])max, which is mean 
metalhcity of the most luminous GCs in the population. The values of these 
parameters in dex are displayed in the upper left of each panel; for reasons 
given in the text, we show results for only two values of ([Fe/H])n,ax- At any 
given luminosity, the GCs follow a unimodal metalhcity distribution of width 
""[Fc/H] = 0-5 dex. Here we use ctlf = 1-4 mag for the GCLF width, appro- 
priate for massive early-type galaxies. We note that several of these panels 
show a strong preference for "blue tilts," i.e., negative blue-peak slopes that 
are significantly larger in absolute value than the red-peak slopes. 



luminosity; thus, these are the populations that would empir- 
ically be described as having a positive MMR. The third row 
of panels in Figure |4] is particularly interesting, as the blue- 
peak slopes are generally larger in absolute value than the 
red-peak slopes. For example, in the three panels of Figure|4] 
with ([Fe/H])„ax = 0.0 and ([Fe/H])o = -0.5, -0.6, -0.7 dex, 
the mean {g-z) blue-peak slopes are {d{g-z)/dz)B = -0.023, 
-0.039, -0.046, respectively, and the red-peak slopes are 
{d{g-z)/dz)R = -0.001, -0.014, -0.024. The scatter in these 
slope values is typically ^ 0.02, so the uncertainties on the 
mean slopes are about ±0.002. If a linear transformation 
were used to convert the slopes from color to metallicity, 
the conclusion would be that the blue GCs have a stronger 
MMR, but this is not the case here, where there is a sin- 
gle MMR and no bimodality in metallicity. For compari- 
son, analysis of a combination ACS Virgo and Fornax Cluster 
Survey (Cote et al. 2004; Jordan et al. 2007a) galaxies gives 
mean (g-z) blue and red peak slopes of -0.035 ± 0.008 and 
-0.016 ± 0.010, respectively, with substantial scatter among 
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Fig. 5. — Same as Fig. |4] but here shown for the slopes with mag of the 
fitted KMM blue and red peaks in the simulated GC (B—I) color distributions. 



different subsamples (Mieske et al. 2010). Thus, the simu- 
lation with ([Fe/H])o, {[FdH])^^^ = -0.6,0.0 dex provides a 
very good match to the observed (g-z) blue and red tilts. 

Even for the top two rows of panels in Figure |4] there is a 
general preference for blue tilts over red ones. More specifi- 
cally, where both blue and red slopes are negative, there is a 
tendency for the blue slope to be more negative, that is, to have 
a more pronounced color-magnitude relation. For example, 
the simulation with ([Fe/H])o, ([Fe/H])max = -0.9,-0.3 (sec- 
ond row, middle panel) has mean (g-z) blue and red slopes of 
-0.025 ±0.002 and -0.016 ±0.002, respectively. We choose 
this as a lower metallicity comparison because it has the same 
factor-of-four range in the MMR limits and also favors the 
blue tilt, though to a lesser degree. 

Figure |5] presents analogous results from the same set of 
simulations but for the slopes of the blue and red peaks of 
the (B-I) color distributions with magnitude. As described 
above, the transformation from metallicity to (B-I) colors is 
based on the 13 Gyr Teramo SPoT models, whereas the (g-z) 
colors are from the P06 empirical transformation. Overall, 
the results are fairly similar, but the tendency to favor blue 
tilts is somewhat weaker among these simulations using the 
model-based transformations. In this case, the model with 
([Fe/H])o, ([Fe/H])niax = -0.6,0.0 (second row, middle panel) 
has mean blue and red slopes {d(B-I)/dI)B = -0.049 ±0.003 
and {d(B-I)/dI)R = -0.032 ±0.002. The blue slope is very 
similar to the value of -0.05 found by Harris (2009), but that 
study found no significant red slope, although Bailin & Har- 
ris (2009) state that it should present at some level. For the 
model with ([Fe/H])o, ([Fe/H])„ax = -0.9,-0.3, the blue and 
red slopes are -0.021 and -0.026, respectively, so the red tilt 
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Fig. 6. — Example histograms of (g-z) and (B—I) colors for GCs in two 
simulated populations with values of ([Fe/H])o and ([Fe/H])iniix given in the 
upper right of each panel. The top panels are typical of the GC color distribu- 
tions for intermediate luminosity galaxies, while the lower panels resemble 
the distributions for the most massive early-type galaxies (e.g., Peng et al. 
2006). Both of these simulations show a tendency towards blue tilts for the 
(g-z) colors; only the higher metallicity model tends towards blue tilts for the 
(B—I) colors, which are derived from a model, rather than empirical, color- 
metallicity relation. 




-1 
[Fe/H] (dex) 

Fig. 7. — Metallicity distributions (grey histograms) for the simulations 
whose color histograms are displayed in Figure |6] We plot metallicities 
only for model GCs brighter than the GCLF peak, since there is no mass- 
metallicity scaling at fainter luminosities. Including the fainter half of the 
GCLF would shift the peak value lower by ~0.08 dex. For comparison, the 
dotted curves show Gaussian metallicity distributions with mean = ([Fe/H])o, 
""[Fe/H] = ".5, and normalized to have the same total areas as the histograms. 
The dashed lines mark the mean metallicities of the plotted histograms. 
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Fig. 8. — Color-magnitude diagrams for two simulations discussed in the 
text. The values of ([Fe/H])o and ([Fe/H])niax are shown, and the magnitudes 
are given with respect to the GCLF mean. We plot points for 10,000 GCs in 
each simulation. The large circles show the results from the KMM fits, which 
are based on 20 times the number of GCs as are plotted in these panels. The 
dashed blue and red lines are linear fits to the KMM peak results. These 
diagrams are illustrative, but in real populations there may be unmodeled 
correlations of mean GC metallicity and color with galaxy luminosity, GCLF 
width, or other parameters. 

is greater in this case, although neither is very large. Despite 
this, we consider it an interesting lower-metallicity compari- 
son model on account of the empirically-based (g-z) slopes 
discussed above. 

Figure |6] presents example (g-z) and (B-I) histograms for 
the two models we have been discussing. These histograms 
exemplify the range of color distributions typically found in 
intermediate- to high-luminosity elliptical galaxies. The dis- 
tributions are clearly bimodal in all the panels of this figure. 
In the top panels the blue peak dominates, giving an appear- 
ance very similar to the GC color distributions found by P06 
for Virgo galaxies with Mb ~ -19 mag. On the other hand, 
the roughly equal peaks in the histograms of the lower panels 
of Figure|6]are quite similar to those found for the most mas- 
sive early-type galaxies (e.g., P06; Harris et al. 2006). These 
are the galaxies that also show the blue tilt most strongly. It 
is interesting that galaxies with relatively more dominant blue 
peaks would have weaker blue tilts; this would be puzzling 
for a trivial linear conversion between color and metallicity. 
However, as shown by the slope-slope comparisons in Fig- 
ures |4] and |5] the more realistic color-metallicity conversions 
used here reproduce this observational finding. 

In comparison to the color histograms. Figure |7] shows that 
despite the MMR built into the models, the underlying metal- 
licity distributions still appear unimodal. Thus, the bifurca- 
tion in color results purely from the nonlinearity of the color- 
metallicity conversion, not the varying of GC mean metal- 
licity with luminosity. Figure [8] illustrates the z versus (g-z) 
color-magnitude diagrams for these baseline models. For clar- 
ity, the plotted simulations contain 10,000 GCs each. The 
larg er ci rcles show results from KMM fits as described in 
Sec. 12.41 In order to delineate the trends more clearly, the 
KMM fits in Figure |8] are based on 20 times the number of 
GCs as are plotted in the panels. As with the color histograms, 
these color-magnitude diagrams bear a good resemblance to 
those found observationally (e.g., P06, Mieske et al. 2006). 

3.2. Dwarf Galaxies: the Effect of the GCLF 

Observationally, the color-magnitude relation of blue-peak 
GCs is mainly a property of giant ellipticals, and the trend is 
weak or absent in lower luminosity galaxies (e.g. DeGraaff 




g—z blue slope 

Fig. 9. — Same as Fig. |4] but here shown for sets of simulations using 
o"LF = 0.9 mag, typical of the GCLF widths of dwarf galaxies. Although some 
individual realizations may have large "tilts" (slopes), the scatter is large and 
encompasses (0,0) at the < la level. Thus, there is no significant tendency 
for either red or blue tilts. 



et al. 2007; Cantiello et al. 2007). However, Mieske et al. 
(2006, 2010) have detected it in pooled samples of GCs from 
faint early-type galaxies. As we have shown, the simulations 
with lower metallicity, which are more appropriate to low lu- 
minosity ellipticals, have less of a tendency to favor the blue 
tilt than do the higher metallicity simulations. However, the 
likely more important reason for the lack of strong blue tilts 
in dwarf galaxies has to do with their GCLFs. Not only do 
they have many fewer GCs, but their luminosity functions are 
narrower. If the GCLF is too narrow, the GCs will not reach 
high enough luminosities to evince a significant MMR. 

Jordan et al. (2006) found that the width of the GCLF de- 
creases linearly with magnitude. For early-type galaxies with 
Mb > -18, the width ctlf ^ 0.9 mag; this includes the major- 
ity of the ACS Virgo and Fornax Survey galaxies, which are 
actually dominated by dwarf ellipticals. Figure |9] shows re- 
sults for the slopes with magnitude of the peaks in the (g-z) 
color distributions analogous to Figure |4l but for populations 
with a narrower GCLF ctlf = 0.9 mag, appropriate to dwarf 
ellipticals. The scatter in the slope values is roughly a fac- 
tor of two larger for the case with ctlf = 0.9 as compared 
to ctlf = 1.4, and encompasses (0,0) at the < Icr level in all 
cases. Measurements of the slopes in individual realizations 
are generally not significant at more than the 1 a level. Thus, 
the lower metallicities and narrower GCLF can explain why 
the observed GC color-magnitude trends are very weak in low 
luminosity galaxies. 
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(B-I) (B-I) 



Fig. 10. — Similar to Fig. [3] but here based on the color-metallicity predic- 
tions from the 11 Gyr SPoT models. Compared to the 13 Gyr models, there 
are more significant differences in the color distributions derived here from 
the polynomial and spline fits. A more complete set of model metallicities 
would help alleviate this interpolation uncertainty. 



3.3. Age Issues 

So far the discussion has centered on GC color distribu- 
tions from simulations using the empirical (§-z)-metallicity 
relation and a theoretical (B-/)-metallicity relation from the 
13 Gyr SPoT models. As shown previously in Figures |2] 
and |3 both of these produce clearly bimodal color distribu- 
tions from single Gaussians in metallicity. As another exam- 
ple. Figure [TO] shows predicted (B-I) distributions similar to 
those in Figure |3] but now based on the 1 1 Gyr SPoT mod- 
els. The agreement between the distributions derived with the 
polynomial and spline interpolations is worse here than for the 
13 Gyr models, mainly because the polynomial fit is poorer 
At the highest metallicities, the polynomial-derived curves in 
Figure [To] appear almost trimodal, whereas the spline-derived 
curves appear fairly similar to those in Figure |3] except that 
the blue peak is a bit broader. 

To examine the effects of these more irregular color distri- 
butions. Figure [TT] shows results for the slopes of the peaks 
in the (B-I) color distributions analogous to Figure |5] but 
based on the polynomial fits to the color-metallicity relations 
from the 1 1 Gyr SPoT models. The scatter is higher by at 
least ^ 50%, and in some cases several times larger, because 
the distributions are less clearly bimodal, making the KMM 
fits much less robust. As in the case for the 13 Gyr mod- 
els, there is no clear preference for blue tilts in the two rows 
of panels in Figure [TT] with ([Fe/H])n,ax = -0.3 dex. How- 
ever, there is a tendency towards blue tilts in the third row of 
panels, in the sense that the slopes are negative and the blue 
slopes tend to be greater in absolute value than the red slopes. 
In particular, for our fiducial "giant galaxy" model with 
([Fe/H])o, ([Fe/H])n,ax = -0.6,0.0 dex, the mean blue and red 
slopes are (c/(B-/)/t//)B = -0.059 ±0.005 and {d{B-I) / dl) r = 
-0.022 ± 0.003. This indicates an even stronger blue tilt than 
for the 13 Gyr models, though again with much larger scat- 
ter Unfortunately, the behavior of the color-metallicity re- 
lation as a function of age is not well calibrated, and as we 
have seen, the results in this case depend on the interpolation 
method. Of course, the ig-z) results based on the empiri- 
cal color-metallicity transformation do not explicitly assume 
any age for the GCs, and the simulated {g-z) distributions are 
strongly bimodal with a clear tendency for the blue tilt. 

We further note that since the point of inflection of the 
nonlinear color-metallicity relation shifts with age, it is also 
possible to introduce a blue tilt by making the brightest GCs 
have systematically older ages (S.-J. Yoon, private commu- 
nication). However, there is no a priori reason to assume a 
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Fig. 11. — Same as Fig. [5] but here based on color-metallicity relations 
from the 1 1 Gyr models (using the polynomial fits to the color-metallicity 
relation, as in Fig.[5j. 



gradient in age with luminosity, whereas the MMR that we 
have used in our simulations is reasonable based on GC self- 
enrichment. Related to this, it would be very useful to have 
observations of GC color distributions at various lookback 
times to see how the color bimodality changes with time. If 
there is underlying bimodality in metallicity, it should change 
very little. On the other hand, if it is due to the inflection 
in the color-metallicity relation, then the bimodality should 
decrease even at modest lookback times of a few Gyr (e.g., 
compare Figures [3] and [TO]). In practice, the observational re- 
quirements are forbidding, since it involves measuring precise 
GC color distributions at z > 0.2. 

Kalirai et al. (2008) attempted to measure evolution in the 
GC colors of a z=0.089 elliptical galaxy serendipitously ob- 
served in extremely deep //ST imaging of a foreground Galac- 
tic globular cluster. These authors detected two peaks in 
the color distribution and interpreted them as the usual local 
peaks, but shifted substantially redder at z=0.089 (even after 
extinction and k corrections). Large fl«f/-evolution in color at 
a lookback time of just over 1 Gyr would not be easy to under- 
stand in the context of either explanation for color bimodality. 
However, from Fig. 2 of Kalirai et al., it appears that the true 
blue and red peaks are probably where they should be, but 
they are unresolved because of observational error at these 
faint magnitudes. The apparent red peak in the color distri- 
bution is likely contamination from the coolest members of 
the white dwarf sequence in the foreground cluster Thus, this 
study illustrates the difficulty of such observations. However, 
it would be worth trying this experiment with the next gener- 
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ation of astronomical observatories, such as the James Webb 
Space Telescope {JWST) and ground-based extremely large 
telescopes (ELTs). In the meantime, further refinements in 
the model color-metallicity predictions for different ages and 
elemental abundance ratios are greatly needed. 



4. DISCUSSION: BIMODALITY OR NONLINEARITY? 

The nonlinear dependence of GC optical colors on metal- 
licity is found both empirically and from detailed SSP mod- 
eling. Because the observed and predicted relations naturally 
produce bimodal color distributions, this nonlinearity is the 
simplest explanation for the observed universality of bimodal 
color distributions. The simplicity of this explanation is the 
strongest argument in its favor In contrast, universal bi- 
modality in the underlying GC metallicity distributions would 
be difficult to understand in the context of hierarchical forma- 
tion models. However, this argument is purely a theoretical 
bias, and it is important to critically examine the available 
data. Of course, even if the nonlinearity is the principal cause 
of the characteristic double-peaked GC color histograms, this 
would not preclude a predominance of two distinct metallicity 
components in some galaxies. 

For the models explored in the present work, there is no 
bimodality in metallicity, just a subtle increase in the mean 
GC metallicity with mass, as may be expected from self- 
enrichment models. We have done this precisely to test 
whether or not the 'blue tilt" can occur within unimodal 
metallicity distributions. As we have shown, the resulting 
color distributions and color-magnitude relations are simi- 
lar to those found in real galaxies. The correlation between 
luminosities and metallicities in these models is dominated 
by scatter, as shown in Figure [121 which plots the z-band 
magnitude-metallicity relation that would be observed for our 
favored giant elliptical model. The plotted "observed metal- 
licities" in the right panel of Figure [12] are found by using 
the colors that include simulated observational error and then 
inverting the color-metallicity relation. The resemblance of 
this figure to the corresponding one observed by Mieske et al. 
(2006, see their Fig. 12) is striking. 

Very little is actually known about the detailed metallicity 
distributions of GCs in giant ellipticals. In support of their 
hypothesis, Yoon et al. (2006) showed that the histogram of 
Mg/7 values for 150 GCs in M87 from Cohen et al. (1998) 
is reasonably consistent with a Gaussian metallicity distribu- 
tion and the conversions predicted by their SSP models; how- 
ever, this represents only 1 % of the GCs and may not be rep- 
resentative of the full population. In contrast, Strader et al. 
(2007) found that metallicities derived from the Lick index 
measurements for 47 GCs in M49 (Cohen et al. 2003) do ap- 
pear bimodal, in the sense that a double Gaussian provided 
a significantly better fit to their derived metallicities than a 
single Gaussian. A skewed distribution, such as found for 
halo stars in other galaxies (e.g., Harris & Harris 2002; Harris 
et al. 2007) and derived by Yoon et al. (2010) from multi-band 
optical photometry of extragalactic GC systems, would also 
favor two-Gaussian models, even if the enrichment process 
was continuous. Perhaps more importantly, this 47-object 
spectroscopic sample represents only ^ 0.6% of M49's to- 
tal population. It may be that such a small sample of bright 
GCs observed in two spectroscopic masks in two regions of 
M49 shows metallicity bimodality, but that the full population 
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Fig. 12. — Magnitude-metallicity relation for a simulated population of 
10,000 GCs with metallicity range parameters (([Fe/H])o, ([Fe/H]),„ax) = 
(-0.6,0.0) and ctlf = 1.4 mag, our fiducial "giant galaxy" model. Left: true 
z-band magnitudes and metallicities for GCs in the simulation. Right: model 
"observed" magnitudes and the metallicities inferred from the (g—z) colors, 
including simulated observational error Consistent with empirical studies, 
there is no well-defined trend of the observed magnitude with metallicity in 
the overall population; the right panel can be compared to Fig. 12 of Mieske 
et al. (2006), to which it bears a strong resemblance. The dashed fiducial 
lines are the same as those in Mieske et al.: the horizontal line marks the 
approximate magnitude where the "tilf becomes more pronounced, and the 
vertical fine approximates the mean metallicity of the GCs fainter than this. 

would have a more continuous distribution from the blending 
of numerous substructures and metallicity components as ex- 
pected from hierarchical accretion. We note that the metallic- 
ities as published by Cohen et al. for the same sample of GCs 
do not strongly favor bimodality (KMM /9-value = 0.10), even 
though they come from the same data and correlate well with 
those of Strader et al. (2007). This is because subtle changes 
in a few metallicities can have large effects on the interpre- 
tation when the sample size is small. In another study, Puzia 
et al. (2005) presented spectroscopic metallicities, ages, and 
[a/Fe] ratios for 71 GCs in seven early-type galaxies, finding 
evidence for trends of decreasing GC age and [a/Fe] with in- 
creasing metallicity. The metallicity histogram for this sample 
appeared unimodal, but with broad tails (see Fig. 7 of Puzia 
et al. 2005). In short, published spectroscopic data sets are 
inadequate for detailed constraints on the GC metallicity dis- 
tributions in the general population of early-type galaxies. 

The prospects for obtaining large samples are of course 
more favorable with multi-band imaging, and the combina- 
tion of high-quality near-IR and optical photometry can pro- 
vide useful constraints on metallicity (see CB07). For in- 
stance, by combining optical and IRAC 3.6/im-band pho- 
tometry for 146 GCs in the peculiar, post-merger elliptical 
NGC5128 (Cen A), Spitler et al. (2008) found evidence for 
bimodality in metallicity. A similar data set for 75 GCs in the 
Sombrero galaxy was inadequate to draw any conclusions in 
this regard. Several other important, ground-breaking studies 
of the optical-IR colors of GCs in giant ellipticals have been 
carried out (e.g., Kissler-Patig et al. 2002; Puzia et al. 2002; 
Hempel et al. 2007). The color distributions generally appear 
broad, often with non-Gaussian tails, but the samples are still 
relatively small and do not show the same obvious bimodal- 
ity as found in the large samples of purely optical GC colors. 
A claimed exception to this was based on matched WFPC2 
/ and NICMOS H photometry of a pruned subsample of 80 
GCs in M87 (Kundu & Zepf 2007). The positions and pho- 
tometry for these objects were never published, but the his- 
togram showed a gap at intermediate color However, the red 
"peak" comprised an excess of only 4 GCs with respect to 
this gap, and a double Gaussian was preferred over a single 
Gaussian at only the 2.4-(t level. Most importantly, the whole 
sample constituted only 0.5% of the full population of M87 
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GCs. Representative samples of optical-lR GC colors in a 
statistical number of elliptical galaxies are needed for general 
constraints on their metallicity distributions. In addition, the 
establishment of empirical conversions from these colors to 
metallicity, and detailed comparison to model predictions as a 
function of age, will be essential. Several groups are actively 
working on this topic, so the data required for these efforts 
should be available in the near future. 

5. CONCLUSIONS 

The combination of a realistic, nonlinear color-metallicity 
dependence with a simple mass-metallicity scaling relation 
motivated by self-enrichment considerations can produce GC 
populations with bimodal color distributions and an apparent 
color-magnitude relation ("blue tilt") for the GCs within the 
blue peak of the color distribution. The underlying metallic- 
ity distributions are unimodal. For certain metaJlicity ranges, 
the color-magnitude relation for GCs within the red peak is 
weak or absent, consistent with observations for most giant 
elliptical galaxies. The models which reproduce these fea- 
tures have nearly equal red and blue peaks, again consistent 
with the data for giant elUpticals. As noted by Yoon et al. 
(2006), this equality occurs because the inflection point in the 
color-magnitude relation is approximately equal to the peak 
metallicity of the GC system. It is when the peak metallic- 
ity begins to move through the inflection point that the "blue 
tilt" becomes most pronounced. As a result, the significance 
of the tilt tends to increase in higher metaUicity populations 
with redder mean colors; thus, it would be greater in more 
massive galaxies and at smaller galactocentric radii, both of 
which are found observationally (Mieske et al. 2006, 2010). 

At the low-metallicity extreme of our simulations, there is 
actually a preference for red tilts over blue ones because the 
peak of the metallicity distribution is weU below the inflec- 
tion point of the color-metallicity relation. But in practice, 
the galaxies that would harbor such metal-poor systems would 
have low luminosities, few GCs, and far too few red GCs to 
show a significant tilt. The narrower GCLFs typical of low- 
luminosity galaxies also greatly reduce the significance of any 
color-magnitude trends in our simulations because of the re- 
duced range in GC luminosity. Further, because we have as- 
sumed a single mass-metallicity relation for all GCs, this sce- 
nario implies that the redder, more metal-rich GC systems of 
giant galaxies will have a higher mean mass than those of 
dwarf elhpticals, again consistent with observations (Jordan 
et al. 2007b). 

A few giant galaxies with weak or absent blue tilts such 
as M49 (Mieske et al. 2006; Strader et al. 2006; but see also 
Lee et al. 2008) can be understood as simply being consis- 
tent with stochastic variations among different model realiza- 
tions. This is because the actual underlying metaUicity dis- 
tribution is broad. If the colors linearly reflected metallic- 
ity, so that the blue peak of the GC color distribution truly 
represented a distinct subpopulation "squeezed" into a nar- 
row metallicity margin that tightly became more metal-rich 
with luminosity, then stochastic variations in such a homoge- 
nous population would be small. Thus, a large galaxy with- 
out a blue tilt would be truly anomalous. However, even 
in our favored high-metallicity model that accurately repro- 
duces the observed mean {g-z} blue and red tilts, at least 15% 
of the simulations lack significant color-magnitude slopes in 
the blue peak. Thus, it is not surprising to find a galaxy like 
M49, and this scenario predicts that others will be found once 
large samples of giant ellipticals with similarly well-measured 



color-magnitude diagrams are available. 

Further observations are needed to test this scenario more 
quantitatively. In particular, we still do not know in detail 
the underlying form of the metaUicity distributions of GCs 
in giant ellipticals, and we have assumed Gaussians only for 
the sake of simplicity, whereas the true distributions may be 
skewed in a manner similar to that of old halo stars (Yoon 
et al. 2010). Accurate characterization of these distributions 
will require still larger spectroscopic samples and improved 
models for converting measured indices to actual metaUici- 
ties. As we have noted, large samples of optical-IR colors in 
multiple elliptical galaxies, and empirical conversions from 
these colors to metallicity, are also needed. The required data 
should be available soon from ongoing observational efforts. 

We have found that the precise shapes of GC color distri- 
butions can depend on the fine detaUs of the color-metalUcity 
conversion. In addition to more extensive data samples, im- 
proving the accuracy of photometrically estimated GC metal- 
licities will require a broader and more densely calculated set 
of models, all with better calibrations against resolved stellar 
populations. The variation of alpha-element enhancement as 
a function of metaUicity needs to be calibrated and taken into 
account for accurate data-model comparisons. Further refine- 
ment in our understanding of changes in the color-metallicity 
relation with age, particularly with regards to the behavior 
of the horizontal branch, is essential. It is also important 
to obtain useful empirical consttaints on trends in GC age 
with metalUcity, although there will likely be large variations 
in any such trends due to the merging histories of individ- 
ual galaxies. A still more challenging task is to measure the 
change in color distributions with age by obtaining high qual- 
ity photometry of GC populations at significant cosmological 
lookback times. 

Finally, we note that this investigation has explored the be- 
havior of the blue tilt phenomenon under variations of only 
two parameters, here labeled ([Fe/H])o and ([Fe/H])max- Al- 
though we consider the results presented here promising, 
more sophisticated models are of course desirable, especially 
once more precise color-metalUcity transformations become 
available. Given the remaining uncertainties and the impor- 
tance of GC metallicity distributions to understanding de- 
tailed galaxy formation histories, we hope that this will con- 
tinue to be an active area of theoretical and observational re- 
search into the coming era of JWST and the planned ELTs. 
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